{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "ca9dc78d",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 32418 (\\N{CJK UNIFIED IDEOGRAPH-7EA2}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 29699 (\\N{CJK UNIFIED IDEOGRAPH-7403}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 24635 (\\N{CJK UNIFIED IDEOGRAPH-603B}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 21644 (\\N{CJK UNIFIED IDEOGRAPH-548C}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 33258 (\\N{CJK UNIFIED IDEOGRAPH-81EA}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 30456 (\\N{CJK UNIFIED IDEOGRAPH-76F8}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 20851 (\\N{CJK UNIFIED IDEOGRAPH-5173}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 20998 (\\N{CJK UNIFIED IDEOGRAPH-5206}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n",
      "/opt/homebrew/Caskroom/miniforge/base/envs/xquant/lib/python3.12/site-packages/IPython/core/pylabtools.py:170: UserWarning: Glyph 26512 (\\N{CJK UNIFIED IDEOGRAPH-6790}) missing from font(s) DejaVu Sans.\n",
      "  fig.canvas.print_figure(bytes_io, **kw)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA08AAAF2CAYAAAC21KNWAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjUsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvWftoOwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAL4FJREFUeJzt3QuYVXW9P/7PcL/UcJG7otxMJUUMAjFLCw6g1lEzE4+GEuHPWx3FG/QPDLU4WscfWTyRlZpPmqYnzUsHJdT8mSiKmWnIEdMAZbiIMAICCvN/vqtmDiPDsFD3zB7m9Xqe9exZl71mrT1r1t7v/f2uzyqpqKioCAAAAGrVpPbZAAAAJMITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADs3yLAQAebzwwgtx2GGHRYsWLWqcv2XLlvjTn/60y2UWLlwYmzZtslwNy/Xt27fG+QAUnvAEwIemoqIihgwZEo899liN8w8//PDcy1iu5uUAqD+67QEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOQgPAEAAOTQLM9CAJDXE088Ee3bt69x3vr163MvY7mdLwdA/SipqKioqKffDQAA0GAUtNveo48+Gl/4wheiR48eUVJSEnffffcun/PII4/EJz7xiWjZsmX069cvbrrpph2WmTlzZvTq1StatWoVQ4cOjfnz5xdoDwAAAOogPG3YsCEOPfTQLOzk8corr8Rxxx0Xn/3sZ+PZZ5+NCy64IL72ta/FAw88ULXM7bffHhMnTozLL788nnnmmWz9o0aNipUrVxZwTwAAgMauzrrtpZanu+66K0444YSdLnPZZZfF/fffH88//3zVtDFjxsTatWtj9uzZ2XhqafrkJz8ZP/rRj7Lxbdu2Rc+ePePrX/96TJo0qQ72BAAAaIyKqmDEvHnzYsSIEdWmpVal1AKVbNmyJRYsWBCTJ0+umt+kSZPsOem5O7N58+ZsqJQC15o1a2KvvfbKQh0AANA4VVRUxFtvvZVdapSyRYMJT2VlZdG1a9dq09J4eXl5vP322/Hmm2/G1q1ba1zmxRdf3Ol6p0+fHtOmTSvYdgMAAA3b0qVLY5999mk44alQUktVuk6q0rp162LffffNXqDS0tI6357/O+d/4qbHX42t23bsMdm0SUmceUSvuPBfPlbn2wUAAI1NeXl5dhnQRz/60V0uW1ThqVu3brFixYpq09J4CjitW7eOpk2bZkNNy6Tn7kyq3JeG90rrrY/wNPaog+IXT6+IJjVcbZZ6EZ5x1EFRWtq2zrcLAAAaq5Icl/MUtNre7ho2bFjMnTu32rQ5c+Zk05MWLVrEoEGDqi2Trl9K45XLNAS9O7WNq08aEE22+/s0LSnJxtP0Xp0EJwAAKDYFbXlKd0NfvHhxtVLkqQR5x44ds25zqTvda6+9FjfffHM2/+yzz86q6F166aXx1a9+NR566KH49a9/nVXgq5S6351xxhkxePDgGDJkSMyYMSMriT5u3LhoSE4e3DMO3rs0jvnBY9n4uCN7xelD9xOcAACgMYanp59+OrtnU6XK645S+Ek3v12+fHksWbKkan7v3r2zoHThhRfGD37wg+yCrZ/97GdZxb1Kp5xySqxatSqmTp2aFZgYOHBgVsb8vUUkGoL99vrfoDTxXz4WbVoUVS9KAACgPu7zVGwXhbVr1y4rHFEf1zxV2rjl3eg/9R83AP7rFaOEJwAAKOJsUFTXPAEAABQr4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAACAH4QkAAKBYwtPMmTOjV69e0apVqxg6dGjMnz9/p8seffTRUVJSssNw3HHHVS1z5pln7jB/9OjRdbErAABAI9Ws0L/g9ttvj4kTJ8asWbOy4DRjxowYNWpULFq0KLp06bLD8r/5zW9iy5YtVeNvvPFGHHrooXHyySdXWy6FpRtvvLFqvGXLlgXeEwAAoDEreMvTtddeGxMmTIhx48ZF//79sxDVpk2buOGGG2pcvmPHjtGtW7eqYc6cOdny7w1PKSxtv1yHDh0KvSsAAEAjVtDwlFqQFixYECNGjPjfX9ikSTY+b968XOv4+c9/HmPGjIm2bdtWm/7II49kLVcHHHBAnHPOOVkLFQAAQIPstrd69erYunVrdO3atdr0NP7iiy/u8vnp2qjnn38+C1Dv7bL3xS9+MXr37h0vv/xyfPOb34xjjjkmC2RNmzbdYT2bN2/Ohkrl5eUfaL8AAIDGp+DXPH0QKTQdcsghMWTIkGrTU0tUpTR/wIAB0bdv36w1avjw4TusZ/r06TFt2rQ62WYAAGDPVNBue506dcpaglasWFFtehpP1ynVZsOGDXHbbbfF+PHjd/l7+vTpk/2uxYsX1zh/8uTJsW7duqph6dKlu7knAABAY1fQ8NSiRYsYNGhQzJ07t2ratm3bsvFhw4bV+tw77rgj62p3+umn7/L3LFu2LLvmqXv37jXOT8UlSktLqw0AAABFVW0vlSn/6U9/Gr/4xS9i4cKFWXGH1KqUqu8lY8eOzVqGauqyd8IJJ8Ree+1Vbfr69evjkksuiSeeeCJeffXVLIgdf/zx0a9fv6wEOgAAQIO85umUU06JVatWxdSpU6OsrCwGDhwYs2fPrioisWTJkqwC3/bSPaAee+yxePDBB3dYX+oG+Nxzz2VhbO3atdGjR48YOXJkXHnlle71BAAAFExJRUVFRTQyqdpeu3btsuuf6rML38Yt70b/qQ9kP//1ilHRpkVR1+8AAIBGnQ0K3m0PAABgTyA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAA5CA8AQAAFEt4mjlzZvTq1StatWoVQ4cOjfnz5+902ZtuuilKSkqqDel526uoqIipU6dG9+7do3Xr1jFixIh46aWX6mBPAACAxqrg4en222+PiRMnxuWXXx7PPPNMHHrooTFq1KhYuXLlTp9TWloay5cvrxr+/ve/V5t/zTXXxHXXXRezZs2KJ598Mtq2bZutc9OmTYXeHQAAoJEqeHi69tprY8KECTFu3Ljo379/FnjatGkTN9xww06fk1qbunXrVjV07dq1WqvTjBkz4lvf+lYcf/zxMWDAgLj55pvj9ddfj7vvvrvQuwMAADRSBQ1PW7ZsiQULFmTd6qp+YZMm2fi8efN2+rz169fHfvvtFz179swC0gsvvFA175VXXomysrJq62zXrl3WHXBn69y8eXOUl5dXGwAAAIomPK1evTq2bt1areUoSeMpANXkgAMOyFqlfvvb38Yvf/nL2LZtWxxxxBGxbNmybH7l83ZnndOnT88CVuWQQhkAAECDrrY3bNiwGDt2bAwcODCOOuqo+M1vfhOdO3eOn/zkJ+97nZMnT45169ZVDUuXLv1QtxkAANjzFTQ8derUKZo2bRorVqyoNj2Np2uZ8mjevHkcdthhsXjx4my88nm7s86WLVtmRSi2HwAAAIomPLVo0SIGDRoUc+fOrZqWuuGl8dTClEfq9veXv/wlK0ue9O7dOwtJ268zXcOUqu7lXScAAMDuahYFlsqUn3HGGTF48OAYMmRIVilvw4YNWfW9JHXR23vvvbPrkpIrrrgiDj/88OjXr1+sXbs2vve972Wlyr/2ta9VVeK74IIL4qqrror9998/C1NTpkyJHj16xAknnFDo3QEAABqpgoenU045JVatWpXd1DYVdEjXMs2ePbuq4MOSJUuyCnyV3nzzzay0eVq2Q4cOWcvV448/npU5r3TppZdmAeyss87KAtaRRx6ZrfO9N9MFAAD4sJRUpBsnNTKpm1+qupeKR9Tn9U8bt7wb/ac+kP381ytGRZsWBc+yAADA+8wGRVdtDwAAoBgJTwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAADkITwAAAMUSnmbOnBm9evWKVq1axdChQ2P+/Pk7XfanP/1pfPrTn44OHTpkw4gRI3ZY/swzz4ySkpJqw+jRo+tgTwAAgMaq4OHp9ttvj4kTJ8bll18ezzzzTBx66KExatSoWLlyZY3LP/LII3HqqafGww8/HPPmzYuePXvGyJEj47XXXqu2XApLy5cvrxp+9atfFXpXAACARqzg4enaa6+NCRMmxLhx46J///4xa9asaNOmTdxwww01Ln/LLbfEueeeGwMHDowDDzwwfvazn8W2bdti7ty51ZZr2bJldOvWrWpIrVQAAAANMjxt2bIlFixYkHW9q/qFTZpk46lVKY+NGzfGO++8Ex07dtyhhapLly5xwAEHxDnnnBNvvPHGh779AAAAlZpFAa1evTq2bt0aXbt2rTY9jb/44ou51nHZZZdFjx49qgWw1GXvi1/8YvTu3Ttefvnl+OY3vxnHHHNMFsiaNm26wzo2b96cDZXKy8s/0H4BAACNT0HD0wf1H//xH3HbbbdlrUyp2ESlMWPGVP18yCGHxIABA6Jv377ZcsOHD99hPdOnT49p06bV2XYDAAB7noJ22+vUqVPWErRixYpq09N4uk6pNt///vez8PTggw9m4ag2ffr0yX7X4sWLa5w/efLkWLduXdWwdOnS97E3AABAY1bQ8NSiRYsYNGhQtWIPlcUfhg0bttPnXXPNNXHllVfG7NmzY/Dgwbv8PcuWLcuueerevXuN81NxidLS0moDAABAUVXbS2XK072bfvGLX8TChQuz4g4bNmzIqu8lY8eOzVqGKl199dUxZcqUrBpfujdUWVlZNqxfvz6bnx4vueSSeOKJJ+LVV1/Ngtjxxx8f/fr1y0qgAwAANMhrnk455ZRYtWpVTJ06NQtBqQR5alGqLCKxZMmSrAJfpR//+MdZlb4vfelL1daT7hP17W9/O+sG+Nxzz2VhbO3atVkxiXQfqNRSlVqYAAAACqGkoqKiIhqZVG2vXbt22fVP9dmFb+OWd6P/1Aeyn/96xaho06Ko63cAAECjzgYF77YHAACwJxCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAchCeAAAAcmiWZyEAAIBXVm+IXz+9NJa9+Xbs06F1fHlwz+jdqW00FsIT0OA09hM3AHXD+011v356aUz6r+eipKQkKioqssef/OHluPqkAXHy4J7RGAhPe4g1G7bEmg2b6+R3VVTUya/hPQr5sjekv+kDL5TF//39/0TJP1+T9DjrDy/HxBEfi5Ef71bfm0cRqijof0/tXnvz7eyYXVG+ObqWtoxRH+8We3doXW/bQ8M/B34YGtv+vl8P/rUsfjD3pR3eby4Y/rH4l/5d63RbXlv7djz4QlmsfGtzdPloy+z9bu/2ret8Gy77r+f+cfxUHkT/fLz0v56LDm1aRI/d2KYWzUqiX5ePRkMjPO0hNmx+N1a9taW+NwMKavm6t7PglM7Vle/9lY/X/v5/Yp8ObaJbu1b1uIXwvx5ZtDKu/39/q/bB644Fy+L/fKZPHPWxLvW9ecAu3m9ScKrp/WbG3P+JfTvW3ftNTeeSO5+p+3PJvX9+vWob3itNv+fPr8epQ/bNvb6WzRtm6QXhCWgwHlm0qtYT98OLVu7WiZvCfOBIf6dV6zdH54+0jKMP6Bzd27VulK9D+rBT0wevnzz6tziga2mdB31/m+q8HjSE95tiOpek/5WKncyr+Of8xqBOIt/MmTOjV69e0apVqxg6dGjMnz+/1uXvuOOOOPDAA7PlDznkkPjd735XbX7qYzl16tTo3r17tG7dOkaMGBEvvfRSgfcCqG9O3MUtfTt60R1/jvueez2e+Nsb2WMa/8P/rIzG+sEravngVbfb42+zPa8HDeX9ppjOJelLhpJatiXNbwwKHp5uv/32mDhxYlx++eXxzDPPxKGHHhqjRo2KlStr/mM//vjjceqpp8b48ePjT3/6U5xwwgnZ8Pzzz1ctc80118R1110Xs2bNiieffDLatm2brXPTpk2F3h0aifRNz6/mL4nrHnope0zj1L9iO3E7Tmr+dnRb+oZ0u8f07WjZusZ1fi6WD17F+Lep7/+bYns9KM7jpFjeb4rpXJJaZytq2ZbPHtA4uiMXvNvetddeGxMmTIhx48Zl4ynw3H///XHDDTfEpEmTdlj+Bz/4QYwePTouueSSbPzKK6+MOXPmxI9+9KPsuanVacaMGfGtb30rjj/++GyZm2++Obp27Rp33313jBkzJve2bdzybjTb8m7Ul/T7a/r5/Xh7y9bY9M7WD2Gr+H8vrYobH3+1Wt/ie597Pb56RO84cv9O9b15jdqwvntlf4uapL/VEX33qrP/g2I6TsrKN2Xb88b6LbHXR1rEp/fvHN1K67ZL2O8Xrqi1i8uchWVx8qDGUYkp6dCmea2vR5pfV8dqMf1tiuH/ppheD4r3OCmW95tiOpekghBfPaJ33PDHV6q2p0nJP7YtTW+/m9uSivl80M+/H5bd2Y6SipRGCmTLli3Rpk2buPPOO7PWo0pnnHFGrF27Nn7729/u8Jx99903a6m64IILqqalVqsUjP785z/H3/72t+jbt2/WKjVw4MCqZY466qhsPIWv99q8eXM2VCovL4+ePXtGzwt+HU1atvmQ9xoAAGgotm3eGEtnfDnWrVsXpaWl9ddtb/Xq1bF169asVWh7abysrKzG56TptS1f+bg765w+fXq0a9euakjBCQAAYHc0imp7kydPzlqz3tvyNP//G77LdNlQLFvzdlZ///3a/M7WOPuWZ7KfZ532iWjZvOmHuHUNZzvS/Rvmv7qmxntglJREDOnVMc4+qm+jek1Sl7Bv3vWXnb4m0088JLrWcRex+lYsx8kdC5bG7OfLsus13it1pRh9cLc6635UbMdJff/fVFpRvike3a5L5Wf271zn/y/F8rcplv+bYnk9iu14zbrK/fHVGrtj1WVX5GI5TopNMZxLPmwtmzeJgT3bRzFI2aD7jCIIT506dYqmTZvGihUrqk1P49261XwzyzS9tuUrH9O0VG1v+2W278a3vZYtW2bDe7Vp0Swb9gStWzSNVh/gZJveTCrd89zrMeKgrvVesjW9eXyQfXo/0omotr7FaX5db1N9vybzXn6j1tfk8ZffaHTlwYvlOHlz4zu1Xryb5tfV8dJrr7bZPUfSBffbX6eQHtP0/fZqG43p/6ZS2u+v1OO+F9Pfplj+b4rl9Sim4zUVY0jXGG3/t6n8UuaGx1+Jg/duV2flsIvlOCk2xXAuKUR4alMkn8Pf3Y3tKGi3vRYtWsSgQYNi7ty5VdO2bduWjQ8bNqzG56Tp2y+fpIIRlcv37t07C1DbL5PSYqq6t7N1suuSrelbuErpm+z6Ktm6fYhL36rXdXUdlWSKu9JPsSiW46RYqkFVSjdrvPbkgfH5AT3i8D57ZY9pvD5uCFvf55JiUwx/m2L5vymW16OYFFM57GI6TqAmBY97qbtcKhAxePDgGDJkSFYpb8OGDVXV98aOHRt77713dl1S8u///u9Z8Yf//M//jOOOOy5uu+22ePrpp+P666/P5peUlGTFJK666qrYf//9szA1ZcqU6NGjR7WiFOx+ydb3fttU1zdfq7yD9vYh7r+fL6vTO2in1rbavpGs65taFoPKD+gVjfy+DsV4nKQPGbVVg6qPDxlp3+u7JbIYziXFqL7/NsXyf1Msr0dtYb+ue38U05dkxXacQJ2Hp1NOOSVWrVqV3dQ2FXRIXetmz55dVfBhyZIl0aTJ/zaAHXHEEXHrrbdmpci/+c1vZgEpVdo7+OCDq5a59NJLswB21llnZVX7jjzyyGyd6aa6NPw7aNd3iEsfrtLvS/teeef59CG0vk7Y9f2mWowf0ItBMRwnPmQU97mE4vy/KTbFEPaL7UsyxwnFrKClyotV6uaXqu7lKUfYUCxdszGWvbn73VLSzefS3dV3dmFm6s7wjc/tH4WWboCX7vC+swvfU5eKYvmGsD7eVCv/PpUX8Nb1N+ipC+fOPqA35m/yi0W6qacPGf/gXEJDC/upm/zO3oNTV8K6+F8ulu2gcWnZvEl8Yt8O0dCyQXFcpcUHli6g7NC2xW4/r3/30pj/tzWxtYbvm5pESTb/kH3aRaHdPG9brfPf2brtQ9uOhvJ9wd/f2BA/3ck36Nc/+rf4woAese9eu3efsve75x/fuzT70HnXn17Lqjru3b51nHjY3rv9+z9sDeRPWXD9e5TG5w4SYpMtW7ftcn56veqF45X3mPPXFdl7bc3vwRF/eW1dfO7A3f/fTjcf3R0Hdf9ofPsLH49p974QJVGSPb/yS7LLv/Dx+OyBnXd7GyicPeW9r0n6RqsBEp72EC2aNcmG3XX64fvFjX98pcZ56eT5lcP3i4+0LPxh0qtT2+x6tprOCGl6ml8X21FM7v9LWa2vyX1/WR6XjT6wzrYnVVtKAxSz3rs4l6T57Vo3r5dtg5qvNaqo9Vqjdm3q5ng944hecdTHOsftTy/NerLs06F1nDK4Z/b+C9RRtT2KX/ogcfVJA7LuLE2blFR7TNPr6qT55cE9d9oilKanE3hjk968antN3k83TdjTOZfQkKSAkoX9GqTpaX5dSu/56Uu5H556WPYoOMGOhCfi5ME946GLjo6zPtMnjhvQI3tM42l6YwtxxaTY3lShIXAuoSER9qHhUTBiDykYsad4dfUGXQb+6ZXVG2L4fz6y0wvfU8BtrK8N7IpzCQ3FHU8vjcv+67nsS7H0kazyMYX9uvwSExqz8t3IBsKT8EQR86YKsOcT9qF+CU+7IDzRkHhTBQAoHKXKYQ9SeQEvAAD1S8EIAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACAHIQnAACA+g5Pa9asidNOOy1KS0ujffv2MX78+Fi/fn2ty3/961+PAw44IFq3bh377rtvfOMb34h169ZVW66kpGSH4bbbbivkrgAAAI1cs0KuPAWn5cuXx5w5c+Kdd96JcePGxVlnnRW33nprjcu//vrr2fD9738/+vfvH3//+9/j7LPPzqbdeeed1Za98cYbY/To0VXjKZwBAAAUSklFRUVFIVa8cOHCLAA99dRTMXjw4Gza7Nmz49hjj41ly5ZFjx49cq3njjvuiNNPPz02bNgQzZr9I+ullqa77rorTjjhhPe1beXl5dGuXbusRSu1igEAAI1T+W5kg4J125s3b17WGlQZnJIRI0ZEkyZN4sknn8y9nsqdqAxOlc4777zo1KlTDBkyJG644YYoUAYEAAAobLe9srKy6NKlS7VpKQB17Ngxm5fH6tWr48orr8y6+m3viiuuiM997nPRpk2bePDBB+Pcc8/NrqVK10fVZPPmzdmwfboEAAAoaHiaNGlSXH311bvssvdBpYBz3HHHZV3/vv3tb1ebN2XKlKqfDzvssKxL3/e+972dhqfp06fHtGnTPvA2AQAAjdduX/O0atWqeOONN2pdpk+fPvHLX/4yLrroonjzzTerpr/77rvRqlWr7DqmE088cafPf+utt2LUqFFZy9J9992XPac2999/f3z+85+PTZs2RcuWLXO1PPXs2dM1TwAA0MiV78Y1T7vd8tS5c+ds2JVhw4bF2rVrY8GCBTFo0KBs2kMPPRTbtm2LoUOH1rrxKTilEHTPPffsMjglzz77bHTo0KHG4JSk6TubBwAAUK/XPB100EFZKfEJEybErFmzslLl559/fowZM6aq0t5rr70Ww4cPj5tvvjkr/JCC08iRI2Pjxo1Zy1Uar7w+KQW2pk2bxr333hsrVqyIww8/PAtWqQz6d7/73bj44osLtSsAAACFvc/TLbfckgWmFJBSlb2TTjoprrvuuqr5KVAtWrQoC0vJM888U1WJr1+/ftXW9corr0SvXr2iefPmMXPmzLjwwguzCntpuWuvvTYLaQAAAA3uPk/FzH2eAACAornPEwAAwJ5EeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAMhBeAIAAKjv8LRmzZo47bTTorS0NNq3bx/jx4+P9evX1/qco48+OkpKSqoNZ599drVllixZEscdd1y0adMmunTpEpdcckm8++67hdwVAACgkWtWyJWn4LR8+fKYM2dOvPPOOzFu3Lg466yz4tZbb631eRMmTIgrrriiajyFpEpbt27NglO3bt3i8ccfz9Y/duzYaN68eXz3u98t5O4AAACNWElFRUVFIVa8cOHC6N+/fzz11FMxePDgbNrs2bPj2GOPjWXLlkWPHj122vI0cODAmDFjRo3z//u//zs+//nPx+uvvx5du3bNps2aNSsuu+yyWLVqVbRo0WKX21ZeXh7t2rWLdevWZa1iAABA41S+G9mgYN325s2bl3XVqwxOyYgRI6JJkybx5JNP1vrcW265JTp16hQHH3xwTJ48OTZu3FhtvYccckhVcEpGjRqV7fQLL7xQ4/o2b96czd9+AAAAKIpue2VlZdn1SNV+WbNm0bFjx2zezvzbv/1b7LffflnL1HPPPZe1KC1atCh+85vfVK13++CUVI7vbL3Tp0+PadOmfQh7BQAANFa7HZ4mTZoUV1999S677L1f6ZqoSqmFqXv37jF8+PB4+eWXo2/fvu9rnan1auLEiVXjqeWpZ8+e73sbAQCAxme3w9NFF10UZ555Zq3L9OnTJyvosHLlymrTU0W8VIEvzctr6NCh2ePixYuz8JSeO3/+/GrLrFixInvc2XpbtmyZDQAAAHUWnjp37pwNuzJs2LBYu3ZtLFiwIAYNGpRNe+ihh2Lbtm1VgSiPZ599NntMLVCV6/3Od76TBbPKboGpml+6uCsVqAAAACiEghWMOOigg2L06NFZ2fHUUvTHP/4xzj///BgzZkxVpb3XXnstDjzwwKqWpNQ178orr8wC16uvvhr33HNPVob8M5/5TAwYMCBbZuTIkVlI+spXvhJ//vOf44EHHohvfetbcd5552ldAgAAGuZNclPVvBSO0jVLqUT5kUceGddff33V/HTvp1QMorKaXioz/vvf/z4LSOl5qYvgSSedFPfee2/Vc5o2bRr33Xdf9phaoU4//fQsYG1/XygAAIAGc5+nYuY+TwAAQNHc5wkAAGBPIjwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADkIDwBAADUd3has2ZNnHbaaVFaWhrt27eP8ePHx/r163e6/KuvvholJSU1DnfccUfVcjXNv+222wq5KwAAQCPXrJArT8Fp+fLlMWfOnHjnnXdi3LhxcdZZZ8Wtt95a4/I9e/bMlt/e9ddfH9/73vfimGOOqTb9xhtvjNGjR1eNp3AGAADQ4MLTwoULY/bs2fHUU0/F4MGDs2k//OEP49hjj43vf//70aNHjx2e07Rp0+jWrVu1aXfddVd8+ctfjo985CPVpqew9N5lAQAAGly3vXnz5mUBpzI4JSNGjIgmTZrEk08+mWsdCxYsiGeffTbr7vde5513XnTq1CmGDBkSN9xwQ1RUVHyo2w8AAFAnLU9lZWXRpUuX6r+sWbPo2LFjNi+Pn//853HQQQfFEUccUW36FVdcEZ/73OeiTZs28eCDD8a5556bXUv1jW98o8b1bN68ORsqlZeXv699AgAAGq/dbnmaNGnSTos6VA4vvvjiB96wt99+O7s2qqZWpylTpsSnPvWpOOyww+Kyyy6LSy+9NLsuamemT58e7dq1qxrStVUAAAAFbXm66KKL4swzz6x1mT59+mTXI61cubLa9HfffTerwJfnWqU777wzNm7cGGPHjt3lskOHDo0rr7wya11q2bLlDvMnT54cEydOrNbyJEABAAAFDU+dO3fOhl0ZNmxYrF27NrtuadCgQdm0hx56KLZt25aFnTxd9v71X/811+9K10V16NChxuCUpOk7mwcAAFCv1zyla5VSKfEJEybErFmzslLl559/fowZM6aq0t5rr70Ww4cPj5tvvjkr/FBp8eLF8eijj8bvfve7HdZ77733xooVK+Lwww+PVq1aZWXQv/vd78bFF19cqF0BAAAo7H2ebrnlliwwpYCUquyddNJJcd1111XNT4Fq0aJFWfe87aXqefvss0+MHDlyh3U2b948Zs6cGRdeeGFWYa9fv35x7bXXZiENAACgUEoqGmGN73TNUyocsW7duigtLa3vzQEAABpANijYfZ4AAAD2JMITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABADsITAABAfYan73znO3HEEUdEmzZton379rmeU1FREVOnTo3u3btH69atY8SIEfHSSy9VW2bNmjVx2mmnRWlpabbe8ePHx/r16wu0FwAAAAUOT1u2bImTTz45zjnnnNzPueaaa+K6666LWbNmxZNPPhlt27aNUaNGxaZNm6qWScHphRdeiDlz5sR9990Xjz76aJx11lkF2gsAAIB/KKlIzT0FdNNNN8UFF1wQa9eurXW5tBk9evSIiy66KC6++OJs2rp166Jr167ZOsaMGRMLFy6M/v37x1NPPRWDBw/Olpk9e3Yce+yxsWzZsuz5eZSXl0e7du2y9acWLAAAoHEq341s0CyKxCuvvBJlZWVZV71KaSeGDh0a8+bNy8JTekxd9SqDU5KWb9KkSdZSdeKJJ9a47s2bN2dDpfTCVL5QAABA41X+z0yQp02paMJTCk5JamnaXhqvnJceu3TpUm1+s2bNomPHjlXL1GT69Okxbdq0Hab37NnzQ9p6AACgIXvrrbeyxpsPLTxNmjQprr766lqXSV3rDjzwwCgmkydPjokTJ1aNb9u2LSs8sddee0VJSUm9J90U4pYuXaoLIUXP8UpD4nilIXG80pCU72HHa2pxSsEpzyVAuxWe0vVIZ555Zq3L9OnTJ96Pbt26ZY8rVqzIqu1VSuMDBw6sWmblypXVnvfuu+9mQajy+TVp2bJlNmwvbwXAupIOvD3h4KNxcLzSkDheaUgcrzQkpXvQ8bqrFqf3FZ46d+6cDYXQu3fvLADNnTu3KiylVJuuZaqs2Dds2LCs8MSCBQti0KBB2bSHHnooa0lK10YBAAA0uFLlS5YsiWeffTZ73Lp1a/ZzGra/J1Pq3nfXXXdlP6fuc6kq31VXXRX33HNP/OUvf4mxY8dmzWcnnHBCtsxBBx0Uo0ePjgkTJsT8+fPjj3/8Y5x//vlZMYm8lfYAAADej4IVjEg3u/3FL35RNX7YYYdljw8//HAcffTR2c+LFi2qqnyXXHrppbFhw4bsvk2phenII4/MSpG3atWqaplbbrklC0zDhw/PquyddNJJ2b2hGqrUnfDyyy/foVshFCPHKw2J45WGxPFKQ9KyER+vBb/PEwAAwJ6gYN32AAAA9iTCEwAAQA7CEwAAQA7CEwAAQA7CUz2bOXNm9OrVK6somO5VlUqwQ7H59re/nd1OYPsh3WoAisGjjz4aX/jCF7JbVqRj8+677642P9VFShVg0w3YW7duHSNGjIiXXnqp3raXxm1Xx+uZZ565w/k23aYF6tr06dPjk5/8ZHz0ox+NLl26ZLcOWrRoUbVlNm3aFOedd17stdde8ZGPfCSrgr1ixYrYkwlP9ej222+PiRMnZqUen3nmmTj00ENj1KhRsXLlyvreNNjBxz/+8Vi+fHnV8Nhjj9X3JkEm3eIinT/Tl1E1ueaaa7JbWsyaNSu78Xrbtm2zc21604diO16TFJa2P9/+6le/qtNthOQPf/hDFoyeeOKJmDNnTrzzzjsxcuTI7BiudOGFF8a9994bd9xxR7b866+/Hl/84hdjT6ZUeT1KLU0p0f/oRz/Kxrdt2xY9e/aMr3/96zFp0qT63jyo1vKUvh1NN7qGYpa+pU83X6+8uXp6i0vf8F900UVx8cUXZ9PS/QW7du0aN910U3aTdSiW47Wy5Snd6/K9LVJQ31atWpW1QKWQ9JnPfCY7l3bu3DluvfXW+NKXvpQt8+KLL8ZBBx0U8+bNi8MPPzz2RFqe6smWLVtiwYIFWfeRSummv2k8HXBQbFI3p/QhtE+fPnHaaafFkiVL6nuTYJdeeeWVKCsrq3aubdeuXfbllXMtxeqRRx7JPqQecMABcc4558Qbb7xR35sEWVhKOnbsmD2mz7GpNWr782vq0r/vvvvu0edX4amerF69OrZu3Zp9+7m9NJ7e6KGYpA+a6Vv62bNnx49//OPsA+mnP/3peOutt+p706BWledT51oaitRl7+abb465c+fG1VdfnX3Lf8wxx2SfGaC+pN5RF1xwQXzqU5+Kgw8+OJuWzqEtWrSI9u3bN6rza7P63gCg+KU37koDBgzIwtR+++0Xv/71r2P8+PH1um0Ae5Ltu5Iecsgh2Tm3b9++WWvU8OHD63XbaLzStU/PP/+86521PNWfTp06RdOmTXeoSJLGu3XrVm/bBXmkb5k+9rGPxeLFi+t7U6BWledT51oaqtRVOn1mcL6lvpx//vlx3333xcMPPxz77LNP1fR0Dk2XoaRr9BrT+VV4qiepmXPQoEFZs/z2TaJpfNiwYfW6bbAr69evj5dffjkr/QzFrHfv3tmb+Pbn2vLy8qzqnnMtDcGyZcuya56cb6lrqeBOCk6pqMlDDz2UnU+3lz7HNm/evNr5NZUyT9dE78nnV9326lEqU37GGWfE4MGDY8iQITFjxoys/OO4cePqe9OgmlSlLN2XJHXVS2VIU3n91HJ66qmn1vemQRbmt/9WPl2TlypDpoua04XLqZ/+VVddFfvvv3/25j9lypSs+Mn2Fc6gGI7XNEybNi27V04K/elLqksvvTT69euXldeHuu6qlyrp/fa3v83u9VR5HVO7du2ye+alx9R1P32eTcduaWlpVjE6Bac9tdJeJpUqp/788Ic/rNh3330rWrRoUTFkyJCKJ554or43CXZwyimnVHTv3j07Tvfee+9sfPHixfW9WZB5+OGH0y03dhjOOOOMbP62bdsqpkyZUtG1a9eKli1bVgwfPrxi0aJF9b3ZNFK1Ha8bN26sGDlyZEXnzp0rmjdvXrHffvtVTJgwoaKsrKy+N5tGqKbjNCIqbrzxxqpl3n777Ypzzz23okOHDhVt2rSpOPHEEyuWL19esSdznycAAIAcXPMEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAACQg/AEAAAQu/b/AxUunijvOziwAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 1000x400 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ADF检验p值: 0.0000\n",
      "推荐窗口长度: 30\n",
      "最优参数组合: OrderedDict({'batch_size': 4, 'bidirectional': True, 'early_stopping_patience': 3, 'hidden_size': 14, 'lr': 0.0005198311220169803, 'num_layers': 3})\n",
      "Using device: mps\n",
      "\n",
      "预测的下一期红球组合: [ 4  9 14 19 23 28]\n",
      "预测的下一期蓝球号码: 8\n",
      "\n",
      "最近几期的真实开奖结果:\n",
      "               date  red1  red2  red3  red4  red5  red6  blue\n",
      "1899  2025-08-14(四)     9    11    12    24    25    26    10\n",
      "1900  2025-08-17(日)    11    13    17    19    23    29    16\n",
      "1901  2025-08-19(二)    15    16    22    23    26    32     4\n",
      "1902  2025-08-21(四)     7     9    11    12    16    29    15\n",
      "1903  2025-08-24(日)     3     5    16    23    26    31    14\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 35757 (\\N{CJK UNIFIED IDEOGRAPH-8BAD}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 32451 (\\N{CJK UNIFIED IDEOGRAPH-7EC3}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 36807 (\\N{CJK UNIFIED IDEOGRAPH-8FC7}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 31243 (\\N{CJK UNIFIED IDEOGRAPH-7A0B}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 20013 (\\N{CJK UNIFIED IDEOGRAPH-4E2D}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 30340 (\\N{CJK UNIFIED IDEOGRAPH-7684}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 25439 (\\N{CJK UNIFIED IDEOGRAPH-635F}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 22833 (\\N{CJK UNIFIED IDEOGRAPH-5931}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 21464 (\\N{CJK UNIFIED IDEOGRAPH-53D8}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n",
      "/var/folders/xj/rt03m8y91wjcylw3wlbcvhzm0000gn/T/ipykernel_67498/2012999428.py:306: UserWarning: Glyph 21270 (\\N{CJK UNIFIED IDEOGRAPH-5316}) missing from font(s) DejaVu Sans.\n",
      "  plt.savefig('model_checkpoints/training_curve.png')\n"
     ]
    }
   ],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import torch\n",
    "import torch.nn as nn\n",
    "from sklearn.preprocessing import MinMaxScaler\n",
    "from sklearn.base import BaseEstimator\n",
    "from skopt import BayesSearchCV\n",
    "from skopt.space import Integer, Real\n",
    "from statsmodels.tsa.stattools import acf, pacf, adfuller\n",
    "from statsmodels.graphics.tsaplots import plot_acf\n",
    "import matplotlib.pyplot as plt\n",
    "import os\n",
    "import json\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "# 创建模型保存目录\n",
    "os.makedirs('model_checkpoints', exist_ok=True)\n",
    "\n",
    "# 1. 数据加载和预处理\n",
    "def load_data():\n",
    "    df = pd.read_csv('./lottery_data.csv')\n",
    "    # 反转数据顺序，使最早的日期在前\n",
    "    df = df.iloc[::-1].reset_index(drop=True)\n",
    "    \n",
    "    red_balls = df[['red1','red2','red3','red4','red5','red6']].values\n",
    "    blue_balls = df['blue'].values.reshape(-1, 1)\n",
    "    combined_data = np.hstack([red_balls, blue_balls])\n",
    "    scaler = MinMaxScaler()\n",
    "    scaled_data = scaler.fit_transform(combined_data)\n",
    "    return scaled_data, scaler, df\n",
    "\n",
    "# 2. 动态窗口长度选择\n",
    "def find_optimal_window(data, max_lag=216, min_lag=30, threshold=0.1):\n",
    "    series = np.sum(data[:, :-1], axis=1)\n",
    "    \n",
    "    # 绘制自相关图\n",
    "    fig, ax = plt.subplots(figsize=(10, 4))\n",
    "    plot_acf(series, lags=20, ax=ax)\n",
    "    plt.title(\"红球总和自相关分析\")\n",
    "    plt.show()\n",
    "    \n",
    "    # 使用偏自相关（PACF）\n",
    "    pacf_values = pacf(series, nlags=max_lag)\n",
    "    \n",
    "    # 检查数据平稳性\n",
    "    adf_result = adfuller(series)\n",
    "    print(f\"ADF检验p值: {adf_result[1]:.4f}\")\n",
    "    if adf_result[1] > 0.05:  # 非平稳数据\n",
    "        diff_series = np.diff(series)  # 一阶差分\n",
    "        pacf_values = pacf(diff_series, nlags=max_lag-1)\n",
    "    \n",
    "    # 改进阈值判断逻辑\n",
    "    for lag, value in enumerate(pacf_values[1:], start=1):\n",
    "        if abs(value) < threshold:\n",
    "            return max(min_lag, lag)\n",
    "    \n",
    "    # 自适应逻辑\n",
    "    if len(data) < 200:\n",
    "        return max(10, len(data)//10)\n",
    "        \n",
    "    return max(min_lag, max_lag//2)\n",
    "\n",
    "# 3. 构建LSTM模型\n",
    "class LotteryPredictor(nn.Module):\n",
    "    def __init__(self, input_size=6, hidden_size=20, output_size=7, num_layers=2, bidirectional=True):\n",
    "        super().__init__()\n",
    "        self.bidirectional = bidirectional\n",
    "        self.num_directions = 2 if bidirectional else 1\n",
    "        \n",
    "        self.lstm = nn.LSTM(\n",
    "            input_size=input_size,\n",
    "            hidden_size=hidden_size,\n",
    "            num_layers=num_layers,\n",
    "            bidirectional=bidirectional,\n",
    "            batch_first=False\n",
    "        )\n",
    "        self.linear = nn.Linear(hidden_size * self.num_directions, output_size)\n",
    "        \n",
    "    def forward(self, x):\n",
    "        if x.shape[0] < 1:\n",
    "            raise ValueError(f\"输入序列长度不足: {x.shape[0]}\")\n",
    "        out, _ = self.lstm(x)\n",
    "        # 通过线性层输出最终预测结果\n",
    "        out = self.linear(out)\n",
    "        return out\n",
    "\n",
    "# 4. 模型包装器（包含早停机制）\n",
    "class ModelWrapper(torch.nn.Module, BaseEstimator):\n",
    "    def __init__(self, \n",
    "                 hidden_size=20, \n",
    "                 num_layers=2, \n",
    "                 bidirectional=True, \n",
    "                 output_size=7, \n",
    "                 lr=1e-3,\n",
    "                 batch_size=2,\n",
    "                 early_stopping_patience=3):\n",
    "        super().__init__()\n",
    "        self.hidden_size = hidden_size\n",
    "        self.num_layers = num_layers\n",
    "        self.bidirectional = bidirectional\n",
    "        self.output_size = output_size\n",
    "        self.lr = lr\n",
    "        self.batch_size = batch_size\n",
    "        self.early_stopping_patience = early_stopping_patience\n",
    "        self.device = torch.device(\"cpu\")\n",
    "        self.best_loss = float('inf')\n",
    "        self.no_improvement = 0\n",
    "        self.train_losses = []\n",
    "        self.val_losses = []\n",
    "            \n",
    "    def fit(self, X, y=None):\n",
    "        # 数据预处理\n",
    "        data_tensor = torch.FloatTensor(X).permute(1, 0, 2).to(self.device)  # (seq_len, batch_size, input_size)\n",
    "        target_tensor = torch.FloatTensor(y).to(self.device)  # (batch_size, output_size)\n",
    "    \n",
    "        # 确保至少有1个训练样本\n",
    "        batch_size = data_tensor.shape[1]\n",
    "        split_idx = max(int(batch_size * 0.8), 1)\n",
    "    \n",
    "        train_data, val_data = data_tensor[:, :split_idx, :], data_tensor[:, split_idx:, :]\n",
    "        train_target, val_target = target_tensor[:split_idx, :], target_tensor[split_idx:, :]\n",
    "    \n",
    "        # 初始化模型\n",
    "        self.model = LotteryPredictor(\n",
    "            input_size=X.shape[2],\n",
    "            hidden_size=self.hidden_size,\n",
    "            num_layers=self.num_layers,\n",
    "            bidirectional=self.bidirectional,\n",
    "            output_size=self.output_size\n",
    "        ).to(self.device)\n",
    "    \n",
    "        criterion = nn.MSELoss()\n",
    "        optimizer = torch.optim.Adam(self.model.parameters(), lr=self.lr)\n",
    "    \n",
    "        # 训练循环\n",
    "        for epoch in range(100):  # 最大训练轮次\n",
    "            # 训练阶段\n",
    "            self.model.train()\n",
    "            optimizer.zero_grad()\n",
    "            predictions = self.model(train_data)\n",
    "            \n",
    "            # 确保predictions不为空\n",
    "            if predictions.shape[0] == 0:\n",
    "                raise ValueError(\"模型输出为空，检查输入数据和模型结构\")\n",
    "                \n",
    "            # 使用最后一个时间步的预测结果\n",
    "            loss = criterion(predictions[-1, :, :], train_target)\n",
    "            loss.backward()\n",
    "            optimizer.step()\n",
    "    \n",
    "            # 验证阶段\n",
    "            self.model.eval()\n",
    "            with torch.no_grad():\n",
    "                val_pred = self.model(val_data)\n",
    "                val_loss = criterion(val_pred[-1, :, :], val_target).item()\n",
    "    \n",
    "            # 记录损失\n",
    "            self.train_losses.append(loss.item())\n",
    "            self.val_losses.append(val_loss)\n",
    "    \n",
    "            # 早停机制\n",
    "            if val_loss < self.best_loss:\n",
    "                self.best_loss = val_loss\n",
    "                self.best_model_weights = self.model.state_dict()\n",
    "                self.no_improvement = 0\n",
    "            else:\n",
    "                self.no_improvement += 1\n",
    "                \n",
    "            if self.no_improvement >= self.early_stopping_patience:\n",
    "                break\n",
    "                \n",
    "        # 加载最佳权重\n",
    "        self.model.load_state_dict(self.best_model_weights)\n",
    "        return self\n",
    "            \n",
    "    def predict(self, X):\n",
    "        self.model.eval()\n",
    "        with torch.no_grad():\n",
    "            # 确保输入数据有正确的维度\n",
    "            if X.ndim == 2:\n",
    "                # 添加序列长度维度，使其变为3D张量\n",
    "                X = np.expand_dims(X, axis=0)\n",
    "            X_tensor = torch.FloatTensor(X).permute(1, 0, 2).to(self.device)\n",
    "            predictions = self.model(X_tensor)[-1, :, :].cpu().numpy()\n",
    "        return predictions\n",
    "    \n",
    "    def get_params(self, deep=True):\n",
    "        return {\n",
    "            'hidden_size': self.hidden_size,\n",
    "            'num_layers': self.num_layers,\n",
    "            'bidirectional': self.bidirectional,\n",
    "            'output_size': self.output_size,\n",
    "            'lr': self.lr,\n",
    "            'batch_size': self.batch_size,\n",
    "            'early_stopping_patience': self.early_stopping_patience\n",
    "        }\n",
    "    \n",
    "    def set_params(self, **params):\n",
    "        for key, value in params.items():\n",
    "            setattr(self, key, value)\n",
    "        return self\n",
    "\n",
    "# 5. 超参数优化模块\n",
    "def optimize_params(X, y):\n",
    "    search_space = {\n",
    "        'hidden_size': Integer(10, 50),\n",
    "        'num_layers': Integer(1, 3),\n",
    "        'bidirectional': [True, False],\n",
    "        'lr': Real(1e-4, 1e-1, 'log-uniform'),\n",
    "        'batch_size': Integer(1, 4),\n",
    "        'early_stopping_patience': Integer(2, 5)\n",
    "    }\n",
    "    \n",
    "    opt = BayesSearchCV(\n",
    "        estimator=ModelWrapper(output_size=7),\n",
    "        search_spaces=search_space,\n",
    "        n_iter=30,\n",
    "        cv=2,\n",
    "        scoring='neg_mean_squared_error',\n",
    "        verbose=0,\n",
    "        n_jobs=1\n",
    "    )\n",
    "    \n",
    "    opt.fit(X, y)\n",
    "    return opt.best_params_, opt.cv_results_\n",
    "\n",
    "# 6. 模型训练函数\n",
    "def train_model(optimal_lag, best_params):\n",
    "    device = torch.device(\"mps\" if torch.backends.mps.is_available() else \"cpu\")\n",
    "    print(f\"Using device: {device}\")\n",
    "    \n",
    "    # 加载数据\n",
    "    scaled_data, scaler, df = load_data()\n",
    "    X, y = prepare_sequences(scaled_data, sequence_length=optimal_lag)\n",
    "    \n",
    "    # 初始化并训练模型\n",
    "    model = ModelWrapper(output_size=7, **best_params)\n",
    "    model.fit(X, y)\n",
    "    \n",
    "    # 保存最佳模型\n",
    "    torch.save({\n",
    "        'model_state_dict': model.best_model_weights,\n",
    "        'params': best_params,\n",
    "        'scaler': scaler,\n",
    "        'optimal_lag': optimal_lag,\n",
    "        'train_losses': model.train_losses,\n",
    "        'val_losses': model.val_losses\n",
    "    }, 'model_checkpoints/best_model.pth')\n",
    "    \n",
    "    return model, scaler, scaled_data, df\n",
    "\n",
    "# 7. 准备训练数据\n",
    "def prepare_sequences(data, sequence_length=10):\n",
    "    if sequence_length < 1:\n",
    "        raise ValueError(\"sequence_length必须≥1\")\n",
    "    windows = []\n",
    "    labels = []\n",
    "    for i in range(len(data)-sequence_length):\n",
    "        windows.append(data[i:i+sequence_length, :-1])\n",
    "        labels.append(data[i+sequence_length])\n",
    "    return np.array(windows), np.array(labels)\n",
    "\n",
    "# 8. 预测函数\n",
    "def predict_next_combination(model, scaler, data, sequence_length):\n",
    "    # 使用与训练时相同的序列长度\n",
    "    last_sequence = data[-sequence_length:]\n",
    "    last_red_balls = last_sequence[:, :-1]\n",
    "    \n",
    "    with torch.no_grad():\n",
    "        prediction = model.predict(last_red_balls)\n",
    "    \n",
    "    predicted_values = scaler.inverse_transform(prediction.reshape(1, -1)).astype(int)\n",
    "    red_balls = np.clip(predicted_values[0][:6], 1, 33)\n",
    "    blue_ball = np.clip(predicted_values[0][6], 1, 16)\n",
    "    \n",
    "    return red_balls, blue_ball\n",
    "\n",
    "# 9. 执行流程\n",
    "# 1. 加载数据\n",
    "scaled_data, _, df = load_data()\n",
    "\n",
    "# 2. 自动选择窗口长度\n",
    "optimal_lag = find_optimal_window(scaled_data, max_lag=216)\n",
    "print(f\"推荐窗口长度: {optimal_lag}\")\n",
    "\n",
    "# 3. 准备数据\n",
    "X, y = prepare_sequences(scaled_data, sequence_length=optimal_lag)\n",
    "\n",
    "# 4. 超参数优化\n",
    "best_params, search_results = optimize_params(X, y)\n",
    "print(f\"最优参数组合: {best_params}\")\n",
    "\n",
    "# 5. 使用最优参数训练模型\n",
    "model, scaler, data, df = train_model(optimal_lag, best_params)\n",
    "\n",
    "# 6. 可视化训练过程\n",
    "def plot_training_curve(train_losses, val_losses):\n",
    "    plt.figure(figsize=(10, 5))\n",
    "    plt.plot(train_losses, label='Train Loss')\n",
    "    plt.plot(val_losses, label='Validation Loss')\n",
    "    plt.xlabel('Epoch')\n",
    "    plt.ylabel('Loss')\n",
    "    plt.legend()\n",
    "    plt.title('训练过程中的损失变化')\n",
    "    plt.savefig('model_checkpoints/training_curve.png')\n",
    "    plt.close()\n",
    "\n",
    "# 7. 保存搜索结果\n",
    "def save_training_artifacts(model, scaler, best_params, train_losses, val_losses, optimal_lag):\n",
    "    # 保存模型\n",
    "    torch.save({\n",
    "        'model_state_dict': model.best_model_weights,\n",
    "        'params': best_params,\n",
    "        'scaler': scaler,\n",
    "        'optimal_lag': optimal_lag,\n",
    "        'train_losses': train_losses,\n",
    "        'val_losses': val_losses\n",
    "    }, 'model_checkpoints/best_model.pth')\n",
    "    \n",
    "    # 保存元数据\n",
    "    metadata = {\n",
    "        'best_params': best_params,\n",
    "        'optimal_lag': optimal_lag,\n",
    "        'train_loss': train_losses[-1],\n",
    "        'val_loss': val_losses[-1],\n",
    "        'best_val_loss': min(val_losses)\n",
    "    }\n",
    "    \n",
    "    with open('model_checkpoints/metadata.json', 'w') as f:\n",
    "        json.dump(metadata, f)\n",
    "\n",
    "# 8. 执行预测\n",
    "red_combination, blue_combination = predict_next_combination(model, scaler, data, optimal_lag)\n",
    "print(f\"\\n预测的下一期红球组合: {np.sort(red_combination)}\")\n",
    "print(f\"预测的下一期蓝球号码: {blue_combination}\")\n",
    "\n",
    "# 9. 显示最近几期的真实开奖结果\n",
    "print(\"\\n最近几期的真实开奖结果:\")\n",
    "print(df.tail(5)[['date', 'red1', 'red2', 'red3', 'red4', 'red5', 'red6', 'blue']])\n",
    "\n",
    "# 10. 保存训练结果\n",
    "save_training_artifacts(model, scaler, best_params, model.train_losses, model.val_losses, optimal_lag)\n",
    "\n",
    "# 11. 可视化训练曲线\n",
    "plot_training_curve(model.train_losses, model.val_losses)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "xquant",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
